* hellwig-create-Bs

*
*
* This file create's the 1000 simulated beta coefficients based off our 3 bootstrapping approaches (delta method is done in hellwig-makeplots.do)
* All created datasets in this file will be used in the next .do file, hellwig-makeplots.do
* ----------------------------------------------------------------------------
cd "~/Dropbox/Benton-Jordan-Philips/final-JOP-replication"
clear
clear matrix
clear mata
set maxvar 80000
use "data/hellwig/ukdata1.dta", clear
keep timevar fdsupp84 fdunem33 fdinfl21 fdpexp17 polevent falkland el79 el83 el87 el92 el97 trade wilscall thatcher blair // keep only used vars
compress
set seed 03858420


/* All approaches below save the parameter estimates in the following datasets:
	--betas-parametric.dta
	--betas-me.dta
	--betas-resid.dta
	
	Here's the coding of each beta:
		b1 = fdsupp84:fdunem33
		b2 = fdsupp84:fdinfl21
		b3 = fdsupp84:fdpexp17
		b4 = fdsupp84:polevent
		b5 = fdsupp84:falkland
		b6 = fdsupp84:el79
		b7 = fdsupp84:el83
		b8 = fdsupp84:el87
		b9 = fdsupp84:el92
		b10 = fdsupp84:el97
		b11 = fdsupp84:trade
		b12 = fdsupp84:wilscall
		b13 = fdsupp84:thatcher 
		b14 = fdsupp84:blair 
		b15 = fdsupp84:_cons
		b16 = HET:trade
		b17 = HET:wilscall
		b18 = HET:thatcher
		b19 = HET:blair
		b20 = HET:_cons
		b21 = HET:L.arch
		b22 = HET:L.garch
*/





* ----------------------------------------------------------------------------
*	APPROACH 1: PARAMETRIC BOOTSTRAP
* ----------------------------------------------------------------------------
arch fdsupp84 fdunem33 fdinfl21 fdpexp17 polevent falkland el79 el83 el87 el92 el97 trade wilscall thatcher blair, robust nolog arch(1) garch(1) het(trade wilscall thatcher blair)
mat B = e(b)
mat VCV = e(V)
loc names: colnames VCV
di "`names'"
mat list B
clear
corr2data b1 b2 b3 b4 b5 b6 b7 b8 b9 b10 b11 b12 b13 b14 b15 b16 b17 b18 b19 b20 b21 b22, n(1000) means(B) cov(VCV)
gen bootcount = _n
compress
save "data/hellwig/betas-parametric.dta", replace
* ----------------------------------------------------------------------------



* ----------------------------------------------------------------------------
*	APPROACH 2: MAXIMUM ENTROPY BOOTSTRAP
* ----------------------------------------------------------------------------
* NOTE: MAKE SURE TO RUN hellwig-me.R first!
clear
clear matrix
clear mata
set maxvar 80000
use "data/hellwig/ukdata1.dta", clear
keep timevar fdsupp84 fdunem33 fdinfl21 fdpexp17 polevent falkland el79 el83 el87 el92 el97 trade wilscall thatcher blair // keep only used vars
compress
set seed 03858420
gen time = _n // to join
preserve // need to csv -> dta
import delimited "data/hellwig/hellwig-me.csv", encoding(ISO-8859-2) clear 
drop yfirstobslastobs v1 // Note: don't need yfirstobslastobs (that's just fdsupp84)
compress
save "data/hellwig/hellwig-me.dta", replace
restore
joinby time using "data/hellwig/hellwig-me.dta"

* plot synthetic y's for SI:
twoway line meboot_y_1 time, color(gray) lpattern(solid) || ///
	line meboot_y_2 time, color(gray) lpattern(solid) || ///
	line meboot_y_3 time, color(gray) lpattern(solid) || ///
	line meboot_y_4 time, color(gray) lpattern(solid) || ///
	line meboot_y_5 time, color(gray) lpattern(solid) || ///
	line meboot_y_6 time, color(gray) lpattern(solid) || ///
	line meboot_y_7 time, color(gray) lpattern(solid) || ///
	line meboot_y_8 time, color(gray) lpattern(solid) || ///
	line meboot_y_9 time, color(gray) lpattern(solid) || ///
	line meboot_y_10 time, color(gray) lpattern(solid) || ///
	line meboot_y_11 time, color(gray) lpattern(solid) || ///
	line meboot_y_12 time, color(gray) lpattern(solid) || ///
	line meboot_y_13 time, color(gray) lpattern(solid) || ///
	line meboot_y_14 time, color(gray) lpattern(solid) || ///
	line meboot_y_15 time, color(gray) lpattern(solid) || ///
	line fdsupp84 time, color(blue) lpattern(solid) lwidth(medthick) legend(off) lcolor(%30)
	*graph export "Bootstrap-Simulation/fall-2023-RR-submission/figures/hellwig-meboot-syntheticplot.pdf", as(pdf) replace

set obs 1000 // to store all beta-hats
forv i = 1/22 { // empty vectors to store B's. k = 22
	gen b`i' = .
}
* loop through the DVs:
timer on 1
loc place = 1
forv i = 1/2500 { // all me's
	noi di in y "On successful simulation #: `place'"
	if `place' == 1001 {
		* STOP; 1000 beta's created
	}
	else {
		cap arch meboot_y_`i' fdunem33 fdinfl21 fdpexp17 polevent falkland el79 el83 el87 el92 el97 trade wilscall thatcher blair, robust nolog arch(1) garch(1) het(trade wilscall thatcher blair)
		if "`e(converged)'" == "1" {
			forv z = 1/22 {
				mat B = e(b)
				qui replace b`z' = B[1,`z'] in `place' // replace empty beta w/ estimates
			}
			loc place = `place' + 1 // move down 1 beta row
		}
		else {
			* do nothing. GARCH failed
		}
	}
}
timer off 1
timer list 1

keep b1-b22
gen bootcount = _n
compress
save "data/hellwig/betas-me.dta", replace
* ----------------------------------------------------------------------------



* ----------------------------------------------------------------------------
*	APPROACH 3: RESIDUAL BOOTSTRAP, using Enders' suggestion of random draws for starting values
* ----------------------------------------------------------------------------
clear
clear matrix
clear mata
set maxvar 80000
use "data/hellwig/ukdata1.dta", clear
keep timevar fdsupp84 fdunem33 fdinfl21 fdpexp17 polevent falkland el79 el83 el87 el92 el97 trade wilscall thatcher blair // keep only used vars
compress
set seed 03858420
do "scripts-programs/bootr2.ado" // same as bootr but preserves any first-obs missing values in master dataset
do "scripts-programs/bootrone.ado" // draws a single boot of y, e for starting values.

arch fdsupp84 fdunem33 fdinfl21 fdpexp17 polevent falkland el79 el83 el87 el92 el97 trade wilscall thatcher blair, robust nolog arch(1) garch(1) het(trade wilscall thatcher blair)
predict nu, resid // nu_t = \sigma_t\epsilon_t
predict sigma2, variance // sigma^2_t
gen epsilon = nu/sqrt(sigma2) // rescale, nu_t/\sigma_t = \epsilon_t
su epsilon // mean-center
replace epsilon = epsilon-r(mean)
gen _sigma2 = . // for storing sigma2. This gets replaced every bootstrap

* loop across bootstraps
forv b = 1/2500 { //2500
	* Create initial values of y and sigma2 at t=0.
	* draw residuals:
	bootrone, resid(epsilon) // this will be the t=0 error which we need for sigma2
	loc residvalue = r(residvalue) 
	* resid for t=1 error for y
	bootrone, resid(epsilon) 
	loc residvaluey = r(residvalue)
	
	* independently from resid, draw bs at same time for y and sigma2
	bootrone, yvalue(fdsupp84) sigma2value(sigma2)
	
	* init sigma2 at t=1
	replace _sigma2 = _b[ARCH:L.arch]*`residvalue'^2*r(sigma2init) + ///
		_b[ARCH:L.garch]*r(sigma2init) + exp(_b[HET:_cons] + ///
		_b[HET:blair]*blair + _b[HET:thatcher]*thatcher + ///
		_b[HET:wilscall]*wilscall + _b[HET:trade]*trade) in 1
	* and sigma2 at t=2, which requires `residvaluey' (resid observed at t=1)
	replace _sigma2 = _b[ARCH:L.arch]*`residvaluey'^2*l._sigma2 + ///
		_b[ARCH:L.garch]*l._sigma2 + exp(_b[HET:_cons] + ///
		_b[HET:blair]*blair + _b[HET:thatcher]*thatcher + ///
		_b[HET:wilscall]*wilscall + _b[HET:trade]*trade) in 2
	* init y at t=1
	gen y_b`b' = _b[fdunem33]*fdunem33 + _b[fdinfl21]*fdinfl21 + ///
		_b[fdpexp17]*fdpexp17 + _b[polevent]*polevent + ///
		_b[falkland]*falkland + _b[el79]*el79 + _b[el83]*el83 + ///
		_b[el87]*el87 + _b[el92]*el92 + _b[el97]*el97 + _b[trade]*trade + ///
		_b[wilscall]*wilscall + _b[thatcher]*thatcher + _b[blair]*blair + ///
		_b[_cons] + `residvaluey'*sqrt(_sigma2) in 1
	
	* draw bs resids for t=2/T
	cap drop boot_e
	bootr2, old(epsilon) newname(boot_e) minT(2)
	
	* recursively gen _sigma2 and y for 2/T
	replace _sigma2 = _b[ARCH:L.arch]*l.boot_e^2*l._sigma2 + ///
		_b[ARCH:L.garch]*l._sigma2 + exp(_b[HET:_cons] + ///
		_b[HET:blair]*blair + _b[HET:thatcher]*thatcher + ///
		_b[HET:wilscall]*wilscall + _b[HET:trade]*trade) in 3/L
	replace y_b`b' =  _b[fdunem33]*fdunem33 + _b[fdinfl21]*fdinfl21 + ///
		_b[fdpexp17]*fdpexp17 + _b[polevent]*polevent + ///
		_b[falkland]*falkland + _b[el79]*el79 + _b[el83]*el83 + ///
		_b[el87]*el87 + _b[el92]*el92 + _b[el97]*el97 + _b[trade]*trade + ///
		_b[wilscall]*wilscall + _b[thatcher]*thatcher + _b[blair]*blair + ///
		_b[_cons] + boot_e*sqrt(_sigma2) in 2/L
	
}

* plot synthetic y's for SI:
twoway line y_b1 time, color(gray) lpattern(solid) || ///
	line y_b2 time, color(gray) lpattern(solid) || ///
	line y_b3 time, color(gray) lpattern(solid) || ///
	line y_b4 time, color(gray) lpattern(solid) || ///
	line y_b5 time, color(gray) lpattern(solid) || ///
	line y_b6 time, color(gray) lpattern(solid) || ///
	line y_b7 time, color(gray) lpattern(solid) || ///
	line y_b8 time, color(gray) lpattern(solid) || ///
	line y_b9 time, color(gray) lpattern(solid) || ///
	line y_b10 time, color(gray) lpattern(solid) || ///
	line y_b11 time, color(gray) lpattern(solid) || ///
	line y_b12 time, color(gray) lpattern(solid) || ///
	line y_b13 time, color(gray) lpattern(solid) || ///
	line y_b14 time, color(gray) lpattern(solid) || ///
	line y_b15 time, color(gray) lpattern(solid) || ///
	line fdsupp84 time, color(blue) lpattern(solid) lwidth(medthick) legend(off) lcolor(%30)
	*graph export "Bootstrap-Simulation/fall-2023-RR-submission/figures/hellwig-rboot-syntheticplot.pdf", as(pdf) replace

* create containers to hold b's:
set obs 1000 // s.t. we can put all the b's
forv i = 1/22 { // empty vectors to store B's. k = 22
	gen b`i' = .
}
* here's the loop:
timer on 1
loc place = 1
forv i = 1/2500 {
	noi di in y "On simulation #: `place'"
	if `place' >= 1001 {
		* STOP; 1000 betas created
	}
	else {
		cap arch y_b`i' fdunem33 fdinfl21 fdpexp17 polevent falkland el79 el83 el87 el92 el97 trade wilscall thatcher blair, robust nolog arch(1) garch(1) het(trade wilscall thatcher blair)
		if "`e(converged)'" == "1" {
			* post the B's to a vector
			forv z = 1/22 {
				mat B = e(b)
				replace b`z' = B[1,`z'] in `place' // replace empty beta w/ estimates
			}
			loc place = `place' + 1 // move down 1 beta row
		}
		else {
			* do nothing, re-estimate using a new Y
		}
	}
}
timer off 1
timer list 1
keep b1-b22
gen bootcount = _n
compress
save "data/hellwig/betas-resid.dta", replace

* ----------------------------------------------------------------------------
